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Summary 


The  purpose  of  this  grant  was  to  construct  a  stereographic  particle  image  velocimeter 
(SPIV)  for  the  measurement  of  three-dimensional  velocity  vectors  on  planar  slices 
through  fluids  defined  by  laser  Ught  sheets.  The  facility  has  been  constructed  from 
components,  requiring  the  use  of  multiple  vendors.  A  stereographic  particle  image 
velocimeter  using  two  2,048  x  2,048  pixel  video  cameras  has  been  constructed  and 
shown  to  work  .  The  attached  paper  by  Liu,  Adrian,  and  Meinhart  (1996)  describes  the 
system  and  demonstrates  its  applicability  to  turbulence  research.  The  system  was  used  in 
1996  to  make  measurements  in  planes  parallel  to  the  wall  of  a  turbulent  boundary  layer, 
at  the  UIUC  boundary  layer  wind  tunnel  in  order  to  study  the  low  speed  streak  spacing  of 
boundary  layers  and  their  contribution  to  the  Reynolds  shear  stress.  It  was  also  installed 
in  the  National  Diagnostic  Facility  of  the  Illinois  Institute  of  Technology  for  use  in  a 
collaborative  study  of  non-equilibrium  turbulent  boundary  layers.  The  system  is  one-of  - 
a-kind,  and  it  has  been  an  unqualified  success. 

A  cost  extension  of  the  grant  until  January  31,  1997  was  requested  and  approved  to 
permit  completion  of  the  entire  system.  The  computer  system  used  for  the  first,  current 
version  of  the  system  was  not  fast  enough  to  make  feasible  the  analysis  of  large  data 
bases  of  vector  data,  and  the  technical  capabilities  of  the  2k  x  2k  cameras  did  not  allow 
for  the  determination  of  the  direction  of  the  velcoity  vectors.  The  extension  was  needed 
to  give  time  to  specify  and  acquire  a  new,  significantly  faster  computer  system  and  new 
camaeras.  The  new  system  has  cross-correlating  cameras  that  permit  the  determination  of 
the  direction  of  the  flow  without  using  pulse  coding. 

Lastly,  a  new  method  of  correcting  for  distortions  which  seriously  impact  the  accuracy  of 
the  stereo  camera  was  developed.  The  method  is  described  in  the  attached  manuscript  by 
Soloff,  Adrian  and  Liu  (1997),  which  will  appear  in  Measurement  Science  and 
Technology.  This  method  and  the  camera  constructed  to  implement  are  now  functional, 
all  programming  and  optical  assembly  having  been  completed.  It  is  currently  undergoing 
shakedown  and  application  in  a  study  of  mixing  by  rotating  turbine  impellors. 

Expenditures  are  consistent  with  the  originally  proposed  list  of  equipment,  see  Table  1. 
To  achieve  the  speeds  needed  for  large  capacity  data  analysis  it  was  necessary  to  re-direct 
some  of  the  funds  to  be  spent  in  the  final  year  of  the  grant  to  the  acquisition  of  the  new 
high  speed  computer,  instead  of  the  Mac  and  two  printers  listed  toward  the  bottom  of 
Table  1.  Permission  to  put  the  remaining  funding  into  the  new  computer  was  requested 
last  year.  All  items  purchased  under  the  equipment  grant  are  summarized  in  Table  1. 

The  University  of  Illinois  has  contributed  in  excess  of  $41,  670.80  in  matching  funds,  and 
TSI,  Inc.  contributed  $25,000  in  software  (see  attached  Report  of  Gift  Received). 
Additional  moneys  in  excess  of  $33,329.20  were  derived  from  non-Federal  grants  and 
contracts.  The  total  contribution  from  these  non-AFOSR  sources  exceeded  the  $100,000 
level  indicated  in  the  original  proposal. 
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ABSTRACT 

Development  of  a  stereoscopic  particle  image 
velocimeter  for  the  measurement  of  three-dimensional 
vectors  on  a  planar  domain  is  described.  The  camera  is 
based  on  two  high-resolution  video  cameras. 
Experiments  in  a  turbulent  boundary  layer  at  ReQ  =  2525 
demonstrate  its  applicability  to  the  measurement  of 
three-dimensional  turbulent  flow.  In  addition  to  the 
quantitative  value  of  the  out-of-plane  component,  it  is 
found  that  having  the  complete  three-dimensional  vector 
also  signiHcantly  improves  the  qualitative  visualization 
of  the  flow. 

1.  INTRODUCTION 

The  extension  of  particle  image  velocimetry  to 
measurement  of  three-dimensional  velocity  vector  fields 
is  both  desirable  and  achievable.  Hinsch  (1995)  reviews 
many  of  the  techniques  used  for  three-dimensional 
measurements.  In  general,  there  are  two  classes  of 
measurement  techniques:  those  which  measure  three- 
dimensional  velocity  vectors  on  full  three-dimensional 
domains,  i.e.  volumes,  and  those  that  measure  three- 
dimensional  vectors  on  planar  domains,  i.e.  light 
sheets.  Three-dimensional  volumetric  measurements  can 
be  performed  by  tising  holographic  techniques  (Barnhart 
et  al.,  1994,  Meng  and  Hussain,  1995),  by  using 
photogrammetric  particle  tracking  techniques  (Maas  et  al 
1993,  Brodkey  1977,  Nishino  et  al  1991),  or  by 
scanning  laser  light  sheets  rapidly  (Brucker  1995). 
Photogrammetric  methods  use  three  or  four  cameras  to 
view  particles  from  several  different  directions. 
Volumetric  measurements  are  in  general  difficult  and 
require  equipment  that  is  rather  different  than  the  standard 
planar  PIV  equipment  using  a  laser  light  sheet  and  a 
single  camera.  Holographic  systems  require  a  unique  set 
of  techniques  and  apparatus,  and  three-dimensional 
particle  tracking  systems  involve  multiple  cameras  and 
specialized  software,  plus  they  place  demands  on  the 
optical  access  required  to  view  the  flow  field.  However, 
the  most  demanding  aspect  of  volumetric  measurements 
is  the  sheer  number  of  velocity  vectors  that  are  obtained 
from  the  measurements.  A  modest  100x100x100 
velocity  grid  measurement  yields  one  million  vectors, 
and  if  ^ese  measurements  are  to  be  repeated  hundreds  or 
even  thousands  of  times  for  the  purpose  of  statistical 
averaging,  the  amount  of  data  is  overwhelming.  Thus, 
while  volumetric  measurements  can  be  extraordinarily 
valuable  for  studying  instantaneous  flow  field  structures. 


they  are  often  too  rich  in  data  to  be  used  on  a  routine 
basis  with  present  equipment 

An  intermediate  approach  to  three-dimensional 
measurements  with  PTV  is  to  perform  measurements  of  a 
three-dimensional  velocity  field  on  a  planar  domain.  In 
this  way,  the  number  of  vectors  is  the  same  as  in 
ordinary  two-dimensional  planar  PIV,  but  the  third 
component  contributes  significantly  to  the 
experimenter’s  capability  to  visualize  the  flow,  and  the 
third  component  is  valuable  for  the  purposes  of 
quantitative  analysis  of  the  flow.  Further,  full 
measurement  of  the  three-dimensional  vectors  eliminates 
the  perspective  error  that  is  inherent  in  monocular  PIV 
systems  (Adrian  1991,  Prasad  and  Adrian  1993).  This 
error  can  be  quite  significant  if  the  flow  has  significant 
out-of-plane  velocity  component  and  the  angular  field- 
of-view  is  not  smaU. 

Systems  operating  with  this  capability  can  be 
designed  using  various  aspects  of  image  correlation 
(Robinson  and  Rockwell  1993),  or  the  change  in  image 
magnification  (Willert  and  Gharib  1992).  However,  the 
most  common  technique  is  the  well-known  stereoscopic 
method  whereby  one  obtains  a  pair  of  images,  each 
image  viewing  the  particles  in  the  illuminated  light 
sheet  plane  from  a  different  direction.  Displacements  of 
the  particle  images  in  the  two  different  views  differ 
because  of  the  different  viewing  angles,  and 
measurements  of  these  displacements  can  be  used  to 
solve  for  the  full  three-dimensional  displacement  of 
particles.  This  technique  can  be  applied  to  individual 
particles,  as  in  particle  tracking  systems,  or  it  can  be 
applied  to  correlation  technique  for  groups  of  particles, 
llie  correlation  approach  is  desirable  because  three- 
dimensional  velocity  vectors  can  also  be  obtained  on 
uniform  grids,  and  the  problem  of  matching  individual 
pairs  is  solved  automatically  accounted.  The  limitation 
of  this  technique  is  that  to  assure  good  measurements  of 
the  out-of-plane  component  of  velocity,  the  viewing 
angles  betwem  the  stereo  lenses  must  be  substantial,  of 
the  order  of  thirty  degrees,  and  the  time  between 
exposures  of  the  images  must  be  small  enough  that  the 
particles  remain  within  the  thickness  of  the  light  sheet 
These  requirements  have  been  explored  amply  and  shown 
to  be  achievable  and  not  inconvenient  in  practice  (Troy 
1994). 

Stereographic  systems  of  conventional  form  have 
been  used  by  Arroyo  and  Created  (1991),  Prasad  and 
Adrian  (1993),  and  Troy  (1994).  These  experiments  each 
used  two  photographic  recordings,  interrogated  the 
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recordings  to  obtain  displacements  and  then  solved  the 
stereo  equations  to  obtain  three-dimensional  velocity 
vectors.  The  difficulty  with  using  photographic  film  is 
the  problem  of  registering  the  two  pieces  of  film  with 
respect  to  each  other  so  that  measurements  of 
displacements  locations  in  the  Hist  image  can  be 
identified  precisely  with  locations  in  the  second.  This 
problem  must  be  faced  for  each  pair  of  photographic 
exposures  that  are  taken  in  a  stereo-photographic 
system.  The  registration  process  is  more  laborious  and 
time  consuming  than  developing  the  film  itself,  and  it 
places  a  fimdamental  constraint  on  the  utility  of  this 
method. 

To  make  feasible  extensive  quantitative  analysis  of 
three-dimensional  vector  fields,  it  is  necessary  to  make 
the  image  acquisition  analysis  process  easy  enough  and 
fast  enough  to  permit  thousands  of  images  to  be  taken 
for  averaging  purposes.  To  this  end,  the  photographic 
recording  must  be  eliminated,  and  replaced  by 
video  graphic  recording.  Video  cameras  in  a  stereo 
system  can  be  registered  during  the  construction  of  the 
stereo  camera  and  by  this  means  they  can  be  aligned  once 
and  for  all,  thereby  eliminating  all  future  efforts  at 
registering  the  two  images.  The  registration  can  be 
accomplished  both  by  means  of  mechanical  alignment 
during  construction,  and  by  electronic  image  processing 
alignment  during  the  experimental  stage.  The  electronic 
alignment  is  necessary  in  situations  where  an  aberrating 
medium  is  placed  between  the  object  plane  and  the 
cameras,  because  the  registration  that  was  originally 
achieved  in  the  construction  of  the  camera  is  distorted  by 
the  medium.  In  this  case,  a  calibration  procedure  can  be 
used  to  correct  for  aberrations,  provided  they  are  not  so 
severe  as  to  completely  destroy  the  images  of  the 
individual  particles. 

In  this  paper,  the  design  and  construction  of  a 
stereographic  camera  based  on  video  recording  is 
descried,  and  some  of  the  procedures  used  in  the 
aberration  calibration  and  correction  process  are 
explained.  The  system  is  applied  to  the  visualization  of 
three-dimensional  vector  field  structure  in  a  low 
Reynolds  number  turbulent  boundary  layer. 

2.  STEREOSCOPIC  PIV  SYSTEM 

The  stereoscopic  camera  is  shown  schematically  in 
Fig.  1.  The  camera  lenses  Li  and  L 2  lie  in  a  plane 
parallel  to  the  object  plane  defined  by  the  laser  light 
sheet.  They  are  offset  laterally  by  di^i  and 
respectively.  Since  the  camera  arrays  also  lie  in  a  plane 
parallel  to  the  object  plane,  the  magnification  of  each 
lens  is  constant,  equal  to 

MQ  =  dildo.  (1) 

where  di  and  do  are  the  image  and  object  distances, 
respectively.  Within  the  angular  field  indicated  by  the 
solid  lines  each  lens  images  point  objects  with 
resolution  and  distortion  that  are  within  the 
specification  of  the  lens.  In  the  object  plane  the  size  of 
the  region  that  can  be  imaged  with  good  resolution  by 
both  lenses,  the  joint  field',  is  determined  by  the 
intersection  of  the  respective  angular  fields.  The  lateral 
offsets  d[^j  and  di^2  chosen  to  be  the  largest  values 


that  still  permit  the  angular  fields  to  overlap  over  the 
desired  region  in  the  object  plane,  i.e.  to  make  the  joint 
field  as  large  as  possible  while  still  keeping  the  angle 
between  the  lenses  large  enough  to  achieve  small 
measurement  error  for  the  z-component  of  velocity. 
Nominally,  d^j  and  dL2  are  equal,  and  the  angle 
between  the  axes  of  the  two  lenses  is  24.4  degrees. 
According  to  Prasad  and  Adrian  (1993)  the  error  in  the  z- 
component  of  velocity  is  about  twice  the  error  in  the  in¬ 
plane  components  of  velocity  when  the  angle  between 
the  axes  of  the  lenses  exceeds  30  degrees. 

The  field  of  view  of  each  CCD  array  is  simply  the 
area  of  the  array  mapped  onto  the  object  plane.  As  in 
Prasad  and  Adrian  (1993)  and  Troy  (1994),  offsetting  die 
cameras  laterally  in  the  image  plane  overlaps  their 
respective  fields  of  view.  Each  camera  sees  nearly  the 
same  region  in  the  fiow,  thereby  maximizing  the  size  of 
the  joint  field  of  view.  This  d^ign  is  relatively  simple 
to  layout,  but  care  must  be  taken  to  ensure  that  tolerances 
in  the  focal  lengths  of  the  lenses  and  in  the  orientations 
of  the  individual  elements  of  the  system  are  taken  into 
account,  either  mechanically,  or  by  careful  calibration. 

The  joint  field  of  view  increases  by  more  than  a 
factor  of  two  if  the  lenses  are  tilted  towards  each  other  so 
as  to  completely  overlay  their  angular  fields.  This 
method,  called  'angular  displacement',  has  been  discussed 
by  Gauthier  and  RiethmuUer  (1988)  and  more  recently  by 
Prasad  and  Jensen  (1996).  Its  principle  advanuge  is  the 
larger  joint  field  of  view,  and  its  principal  disadvantage 
is  the  more  complicated  alignment  that  is  needed,  which 
also  involves  using  magnification  that  varies  across  the 
field  of  view.  Only  the  simpler  lateral  offset  method  will 
be  considered  here,  although  angular  offset  can  be  easily 
accommodated  by  the  present  camera. 
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Figure  1.  Schematic  of  the  stereoscopic  camera. 
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The  camera  uses  two  Kodak  Megaplus  4.2  CCD 
video  cameras,  each  with  2029  x  2044  pixels  resolution- 
Since  the  camera  is  intended  to  image  a  100-300  mm 
5eld  of  view  in  the  flow,  and  the  CCD  chip  dimension  is 
only  18  mm,  the  magniHcation  of  the  system  is 
substantially  less  than  unity.  M  =  0.23  for  the  present 
system.  Two  Nikon  EL-Nikkor  302.5  mm  F/5.6 
barging  lenses  that  are  optimized  for  normal  operation 
at  large  magniHcation  are  ultilized  to  provide  the 
appropriate  small  magnification  with  little  aberration, 
minimal  distortion  and  large  aperture  (  53.6  mm).  The 
angular  field  of  view  of  this  lens  is  57  degrees.  Images 
are  normally  acquired  at  F/9.  The  video  cameras  and  two 
lenses  are  mounted  on  a  3-D  translation  stage.  The  stage 
provides  the  six  degrees  of  freedom  needed  to  position 
and  focus  the  cameras  and  to  accurately  register  them 
with  respect  to  each  other.  For  the  wind  tunnel 
experiments,  Dq  -  1617  mm  and  di  =  372mm. 

The  rest  of  the  stereoscopic  system  is  shown  in 
Fig.  2.  Illumination  of  the  flow  field  is  provided  by  two 
Continuxim  Lasers  NdiYAG  lasers  with  a  wavelength  of 
532  nm  at  pulse  frequency  of  10  Hz.  Each  laser  delivers 
up  to  200  mJ/pulse  of  energy  with  a  pulse  duration  of  8 
nsec.  The  laser  beam  is  formed  into  a  sheet  of  about  1 
mm  thickness  in  the  test  section.  To  resolve  directional 
ambiguity,  image  shifting  is  normally  needed.  However 
image  shifting  was  not  used  for  the  present 
measurements  because  of  the  strong  mean  flow  in 
comparison  with  the  turbulence  intensity  in  the 
boundary  layer.  If  directional  ambiguity  were  a  problem 
the  cameras  could  be  replaced  by  Ik  x  Ik  cross- 
correlating  camera.  Timing  for  the  sequence  of  laser 
pulses  and  the  image  capture  with  the  cameras  is  provided 
by  a  TSI,  Inc.  Laser  Pulse  synchronizer  box. 
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Figure  2,  Stereoscopic  system. 


3.  IMAGE  ANALYSIS 

The  image  analysis  procedure  is  described  in  Fig.  3. 
The  capture  and  analysis  of  images  is  controlled  by 
INSIGHT,  a  Windows-based  software  package  supplied 
by  TSI,  Inc.  The  software  controls  the  simultaneous 
image  capture  by  the  video  cameras.  Image  analysis  to 
determine  the  velocity  vectors  on  each  image  plane  is 
performed  by  the  autocorrelation  technique.  The  outcome 
of  the  analysis  is  a  pair  of  displacement  vector  fields, 


AX j  and  AX 2  that  are  functions  of  the  camera 
coordinates,  Xi  and  The  vectors  in  these  fields  are 
tested  for  validity  using  procedures  described  by 
Meinhart,  et  al,  (1996),  and  invalid  vectors  are  either 


Figure  3.  Image  analysis  procedure. 

replaced  by  second  or  third  choice  vectors  from  the 
autocorrelation  interrogation,  replaced  by  interpolated 
values,  or  simply  removed.  Typically  the  fraction  of 
invalid  vectors  is  less  than  2%.  Corrections  for  the  mis¬ 
registration  of  the  images,  and  the  variation  of  the 
magnification  factors  are  all  performed  in  a  post¬ 
processing  package.  Effectively,  this  step  maps  the 
displacement  vectors  from  a  grid  in  each  of  the  CCD 
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array  image  planes  onto  a  common  grid  in  the  object 
plane.  This  step  is  necessary  in  general  if  there  is  any 
distorting  me^um  between  the  object  plane  and  the 
camera,  such  as  a  water-air  interface  (Prasad  and  Adrian 
1993).  In  the  present  experiment  the  6,35  mm  thick 
glass  windows  of  the  wind  tunnel  cause  the  images  to  be 
displaced  laterally  by  0,55  mm.  Since  this  shift  is  nearly 
constant,  it  is  most  easily  accounted  for  by  mechanic^ 
alignment  in  the  initial  layout  of  the  cameras.  Li  more 
general  cases  of  nonlinear  distortion  it  is  better  to 
calibrate  the  distortion  and  correct  for  it  in  software.  The 
procedures  for  doing  this  will  be  discussed  elsewhere. 


(a) 


shown  that  the  volume  averaged  displacements  are  given 
by 


Lastly,  the  stereoscopic  equations  are  applied  to 
compute  the  third-component  of  the  velocity  vectors. 
For  the  ideal  camera  system  shown  in  Fig.  1  the  image  of 

a  point  particle  at  x  is  mapped  by  lens  (,  onto  a  location 
Til  ^  X^-plane  on  camera  i  given  by 

-Vl  =  +  ^Ll  )X  +  (2a) 

Uq  Z 

-m = 7^[(^  “  + yy]  (2b) 

Up  Z 

where  carats  denote  unit  vectors.  Three-dimensional 
displacement  of  the  particle  by  an  amount  Ax  results  in 
two-dimensional  displacements  AVji  that  are  found  by 

taking  the  increment  of  (2).  The  result  depends  on  the  z- 
locadon  of  the  particle  at  die  time  of  the  first  exposure. 
However,  PIV  correlation  analysis  gives  an  estimate  of 
the  volume  average  of  the  displacements  of  the  particles 
that  lie  in  the  measurement  volume  defined  by  the 
thickness  of  the  light  sheet  and  the  area  of  the 
interrogation  spot.  By  integrating  the  equation  for  the 
image  displacement  over  such  a  volume  and  by  making 
use  of  the  normally  large  ratio  of  the  object  distance  to 
the  light  sheet  thickness,  Prasad  and  Adrian  (1993)  have 


(b) 


Figure  4.  Image  displacement  Helds  due  to  a  pure  out- 
of-plane  translation,  (a)  Camera  1,  — 7Jj(Xj);  (b) 
Camera2,  (c)  Three-dimensional  displacement 

field,  Ax(X). 
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Given  the  two  two-dimensional  vector 
displacements  AT/j,  AT72  three  of  the  equations  in  (3) 
can  be  solved  for  Ax,  Ay^  Az ,  and  the  fourth  can  be  used 
to  obtain  a  second  estimate  of  Ay,  Alternatively,  the 
two  Ay -equations  can  be  added  together  to  obtain  a 
slightly  more  robust  equation.  The  three-dimensional 
velocity  vector  is,  of  course,  found  by  dividing  three 
dimensional  displacement  by  the  time  between 
expostires. 

Figure  4a  and  4b  show  the  displacement  fields  that 
result  on  cameras  1  and  2,  respectively  when  the  object 
^eld  undergoes  a  pure,  uniform  displacement  in  the  z- 
direction.  Solving  equation  (3)  for  the  three-dimensional 
particle  displacement  yields  the  field  shown  in  Fig.  4c. 

4.  TURBULENT  BOUNDARY  LAYER 
MEASUREMENTS 

The  turbulent  boundary  layer  measurements  were 
performed  in  an  eiffel-type  low -turbulence  boundary 
layer  wind  tunnel  with  a  working  section  914  mm  wide  x 
457  mm  high  x  6.096  m  long.  The  free  stream  turbulence 
intensity  at  the  test  section  inlet  is  less  than  0.2%  for 
free  stream  velocities  less  than  10  m/sec.  The  turbulent 
boimdary  layer  was  produced  on  a  flat  plate  placed  100 
mm  above  the  floor  of  the  test  sectioiL  The  boundary 
layer  was  tripped  by  the  a  4.7  mm  diameter  wire  which 
spaimed  the  boundary  layer  plate  just  downstream  of  the 
leading  edge.  Optical  access  to  the  boundary  layer  was 
provided  from  the  side  by  float  glass  windows,  and  from 
below  by  610  mm  wide  x  2.748  m  long  float-glass 
windows.embedded  in  the  boundary  layer  plate.  The 
stereo  PIV  measurements  presented  here  were  performed 
at  a  free  stream  velocity  of  3.4  m/sec,  which  produced  a 
Reynolds  number  based  on  the  momentum  thickness  0  of 
Re  0=2515  at  the  location  of  the  measurements.  All 
measurements  were  made  with  the  light  sheet  in  the  x-y 
plane  of  the  flow,  where  x  is  streamwise  and  y  is  normal 
to  the  walL 

Figure  5  shows  the  components  of  one 

realization  of  the  flow  field.  A  constant  convection 
velocity  of  3000  mm/s  has  been  subtracted  to  make  the 
fluctuations  more  visible.  The  edge  of  the  turbulent 
boundary  layer  lies  around  50  mm  above  the  wall  in  this 
realization.  As  observed  in  earlier  work,  there  is  an 
internal  layer  close  to  the  wall  that  grows  in  a  manner 
similar  to  that  of  an  ordinary  boundary  layer.  This 
internal  layer,  seen  as  a  region  of  uniformly  low 
momentum  that  extends  from  the  wall  up  to  about  20 
mm,  is  capped  by  a  collection  of  several  intense 
spanwise  vortices. 

The  w-component  of  velocity  is  shown  in  Figure  6 
by  plotting  its  contours  in  gray-level  form  .  It  is 
ixnmediately  clear  that  the  values  of  the  w-component  are 
frequently  as  large  or  larger  than  the  other  components. 
There  is  a  tendency  for  the  w-contours  to  align  along  the 
same  30  degree  inclination  that  is  observed  in  the  j<-v 
field.  There  is  also  a  tendency  for  the  sign  of  the  w- 
component  to  be  positive  on  one  side  of  a  30  degree  line 
and  negative  on  the  other,  indicating  a  large-scale 
rotation  about  an  inclined  axis.  The  strength  of  the  w- 
component  and  its  coherent  organization  shows  that  one 


must  be  very  cautious  in  interpreting  the  two- 
dimensional  fields  in  a  data  set  such  as  Figure  5. 


m/s,  m/s  u-U^=300  mm/s 
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Figure  5.  Streamwise-wall  normal  vector  field  of  the 
turbulent  boundary  layer. 
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Figure  6,  Field  of  the  spanwise  turbulent  velocity 
superimposed  on  the  streamwise-wall  normal  field  from 
Fig.  5. 

The  three  dimensional  vectors  are  shown  in  end 
view,  i.e.  looking  at  the  x-y  plane  along  the  x-axis,  in 
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Figure  7a,  and  in  oblique  view  in  Figure  7b.  The  growth 
of  the  internal  low-momentum  layer  is  very  evident  in 
Figure  7b. 
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U^=3.4  m/s,  U^=3  m/s 

TTM  '  C 


(a)  (b) 

Figure  7.  Three-dimensional  vectors  on  the  x-y  plane, 
(a)  end  view;  (b)  oblique  view. 

5.  SUMMARY 

Stereoscopic  PIV  corrects  for  perspective  error  as 
well  as  providing  the  out-of-plane  component.  It  requires 
a  specie  camera,  but  otherwise  it  uses  essentially  the 
same  hardware  as  a  monocular  PIV.  Videographic 
recording  increases  the  ease  of  use  substantially  by 
eliminating  the  necessity  of  registering  photographic 
images.  The  present  system  has  been  used  to  measure 
thousands  of  frames  in  experiments  that  required  only  a 
few  days  to  set  up.  Although  the  video-based  PTV  system 
has  less  resolution  than  a  photographic  system,  it  is 
encouraging  that  it  is  still  able  to  resolve  many  of  the 
important  features  of  the  turbulent  boundary  layer  that 
have  only  recently  been  discovered  by  photographic 
PIV.  The  results  are  very  similar  to  those  obtained  by 
photographic  PIV.  Thus,  while  the  video  system  clearly 
has  a  role  to  play  in  the  acquisition  of  data  for  analysis 
by  statistical  averaging,  it  also  provides  enough 
resolution  to  be  of  considerable  value  in  visualizing  the 
flow. 
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Abstract 

Optical  distortion  caused  by  inaccurate  optical  alignment,  lens  nonlinearity,  and/or 
refraction  by  optical  windows,  fluid  interfaces,  and  other  optical  elements  of  an  experiment  cause 
inaccuracy  by  introducing  variable  magnification.  Since  fractional  changes  in  the  magnification 
have  a  one-to-one  effect  on  the  accuracy  of  measuring  the  velocity,  it  is  important  to  compensate 
for  such  distortions.  A  general  experimental  calibration  procedure  is  described  which  determines 
the  magnification  matrix  of  a  distorted  imaging  system,  and  an  algorithm  is  presented  to  compute 
accurate  velocity  field  displacements  from  measurements  of  distorted  PIV  images.  These 
procedures  form  a  basis  for  generalized  stereoscopic  PIV  procedures  which  permit  easy 
electronic  registration  of  multiple  cameras  and  accurate  recombination  of  stereoscopic 
displacement  fields  to  obtain  the  three-dimensional  velocity  vector  field. 


1.  Introduction 

The  accuracy  of  velocity  measurement  by  particle  image  velocimetry  (PIV)  or  particle 
tracking  velocimetry  (PTV)  depends  upon  the  accuracy  with  which  image  displacements  can  be 
measured  and  the  accuracy  with  which  image  displacements  can  be  related  to  particle  locations 
and  their  respective  particle  displacements.  (A  third  factor,  the  accuracy  with  which  particle 
displacements  are  equal  to  fluid  displacements,  is  not  considered  here.)  These  factors  depend  in 
turn  on  two  types  of  optical  aberrations;  focusing  aberrations  which  cause  the  particle  images  to 
be  enlarged,  making  it  more  difficult  to  measure  displacements  accurately,  and  image  distortion 
which  creates  a  nonlinear  relationship  between  the  location  of  a  particle  and  the  location  of  its 
image.  Focusing  aberrations  include  spherical  aberration,  astigmatism,  and  coma  (Jenkins  and 
White  1957).  Image  distortion  refers  to  any  agency  that  causes  the  magnification  to  vary  as  a 
function  of  position  in  the  field  of  the  image. 

There  are  three  agencies  that  cause  the  magnification  in  an  experiment  to  vary.  The  first 
is  design  imperfection  in  the  lens.  The  second  is  refraction  caused  by  anything,  including 
windows  and  fluid  interfaces,  that  lies  between  the  object  plane  in  the  illuminating  light  sheet 
and  the  image  plane  (c.f.  Prasad  and  Adrian  1993,  for  example).  The  third  is  misalignment  of  the 
image  recording  plane  with  respect  to  the  object  plane.  The  image  recording  plane  is  the  surface 
of  the  film  or  the  video  array  on  which  the  image  is  actually  recorded,  and  it  need  not  coincide 
with  the  image  plane  in  which  all  images  are  in  ideal  focus.  Even  if  they  do  coincide,  the 
magnification  is  not  constant  unless  the  image  plane  is  parallel  to  the  object  plane,  as  in  lateral- 
offset  stereoscopic  cameras,  such  as  those  described  in  Arroyo  and  Created  (1991)  and  Grant  et 
al  (1995).  This  last  distortion  necessarily  occurs  in  angular-offset  stereoscopic  cameras  in  which 
the  image  plane  is  rotated  with  respect  to  the  object  plane  in  order  to  bring  the  images  into  better 
focus  (Prasad  and  Jensen  1995).  It  also  occurs  in  single  lens  cameras  and  in  lateral-offset 
stereoscopic  cameras  due  to  imperfect  alignment  of  the  recording  plane  with  the  object  plane. 

The  PIV  recording  and  interrogation  processes  are  illustrated  in  Fig.  1  for  the  case  of  zero 
distortion.  The  image  of  any  particle  is  simply  located  by  the  chief  ray  that  extends  from  the 
particle  through  the  center  of  the  lens  to  the  image  plane.  Choosing  coordinate  systems  in  the 
object  and  image  planes  with  origins  Ox  and  OX  on  the  optic  axis  of  the  lens,  the  following 
simple  relationship  results; 

(X,Y)=-Mix,y),  (1) 

where  X  =  (X,  F)  is  the  image  location  and  (x,  y)  is  the  object  location  formed  from  the  projection 
of  X  =  (x,  y,  z)  onto  the  xy-plane.  The  magnification  is 

M=-^,  (2) 

do-z 

in  which  di  is  the  image  distance  and  do  is  the  object  distance. 

The  case  of  nonzero  distortion  is  illustrated  in  Fig.  2.  Here,  a  curved  window  acts  like  a 
weak  lens  to  refract  the  chief  ray,  and  the  image  recording  plane  is  tilted  with  respect  to  the 
object  plane.  Both  of  these  effects  cause  some  blurring  of  the  point  response  image  on  the 
recording  plane,  but  as  noted  above,  that  blurring  may  be  small  enough  to  permit  successful 
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Figure  1.  PIV  recording  and  interrogaiion  processes  Tor  a  single  camera  with  zero  distortion.  X/  is  the 
location  ot  the  center  ot  the  interrogation  spot  on  the  image  plane. 
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PLANE 


Figure  2.  PIV  recording  process  tor  a  single  camera  with  nonzero  distortion.  F(x)  is  the  mapping 
runciion  tor  the  camera  that  relates  positions  in  the  obiect  plane  to  positions  m  the  image  plane. 


measurement  of  the  image  displacements  despite  significant  distortion.  The  reference  frame  X  is 
embedded  in  the  recording  frame,  i.e.  Z  =  0,  and  since  it  is  normally  associated  with  the  pixel 
index  array,  its  orientation  in  this  general  case  with  respect  to  the  object  plane  is  allowed  to  be 
arbitrary.  In  general  the  relationship  between  positions  in  the  three-dimensional  physical  space 
of  the  flow,  X,  and  the  positions  in  the  two-dimensional  image  recording  plane,  X,  can  be 
described  by  a  nonlinear  mapping  function, 

X=F(x),  (3) 

that  replaces  the  simple  linear  relationship  in  Eq.  (1).  Note  that  x  in  this  equation  is  the  full, 
three-dimensional  vector. 

Variable  magnification  introduces  systematic  bias  errors  into  the  measurement  of  fluid 
velocity  and  the  location  of  the  velocity  vector.  Being  systematic,  the  distortion  bias  can  be 
corrected  by  an  appropriate  technique  based  either  on  a  mathematical  model  of  the  distortion,  as 
in  Prasad  and  Adrian  (1993),  or  on  an  experimental  calibration,  as  in  the  calibration  for  particle 
tracking  velocimetry  described  by  Kent  et  al.  (1993).  In  the  work  presented  here,  a  calibration  is 
used  to  generate  an  approximate  mapping  function.  Distortion  can  be  compensated  even  when 
the  focusing  aberrations  cause  some  blurring  of  the  images,  so  long  as  the  blurring  is  not  so 
severe  as  to  make  the  measurement  of  image  locations  impossible  or  unacceptably  inaccurate. 

The  information  found  by  calibration  can  be  used  to  do  more  than  correct  for  distortion. 
Most  PrV  interrogation  analyses  determine  the  velocity  vectors  on  a  uniform  grid  fixed  in  the 
pixel  array  of  the  video  camera,  and  the  user  has  no  ability  to  request  data  on  a  non-uniform  grid 
in  the  fluid  that  may  be  more  appropriate  for  the  flow  under  consideration.  A  side  benefit  of 
using  the  mapping  function  algorithm  described  here  is  that  the  measured  vectors  can  be 
interpolated  onto  any  grid  within  the  domain  of  the  mapping  function. 

The  use  of  a  general  mapping  function  found  by  calibration  is  very  helpful  in  the  physical 
construction  of  a  stereoscopic  PIV  camera  because  it  eliminates  the  need  to  mechanically  register 
the  pixel  arrays  of  multiple  cameras,  replacing  this  process  by  a  simple  mathematical  registration. 
In  addition,  it  will  be  shown  how  the  mathematics  of  the  mapping  function  provides  a  natural 
framework  in  which  to  combine  the  images  of  two  or  more  cameras  to  obtain  the  three- 
dimensional  vector  field. 

2.  Displacement  Equations 

2.1  Single  Camera 

Consider  first  the  displacement  of  one  particle  as  recorded  by  a  single  camera.  Following 
this  analysis  the  results  will  be  generalized  to  the  displacement  of  groups  of  particles,  as 
measured  by  PIV  correlation  techniques.  In  the  next  section,  the  results  will  be  applied  to 
stereoscopic  cameras. 

Between  exposures  at  times  t  and  r  +  At  a  particle  at  x  moves  to  x  +  Ax.  The 
displacement  of  its  image  is 

AX=F(x  + Ax) -F(x).  (4) 

To  first  order  the  displacement  is  approximated  by 
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AX  =  VF(x)Ax, 


(5) 


where  the  gradient, 


=  i  =  l,2and;  =  l,2,3, 

OX; 


(6) 


is  clearly  a  generalization  of  the  usual  concept  of  magnification.  We  shall  call  it  the 
magnification  matrix.  As  an  example,  and  for  later  reference,  the  magnification  matrix  for 
perfectly  linear  imaging  is,  from  Eq.  (1), 


Fi.M)  = 


-d^ 


d„-z 


1  0 
0  1 


V 


(7) 


The  terms  in  the  last  column,  Fi,3  and  F2,3,  represent  the  effect  of  perspective. 

In  PTV  the  displacement  of  the  image  of  an  individual  particle,  AX,  is  determined  by 
image  analysis,  and  one  must  then  relate  that  displacement  to  displacement  in  the  physical  space. 
Ax.  Eq.  (4),  or  more  commonly  the  approximate  Eq.  (5),  can  provide  the  necessary  relationship 
when  distortions  are  present,  but  there  are  two  problems  in  executing  the  last  step.  First,  since  X 
is  two-dimensional  and  x  is  three-dimensional,  there  is  not  a  unique  relationship  between  X  and 
X.  The  particle  could  lie  anywhere  along  the  chief  ray  that  extends  from  the  image  point  through 
the  center  of  the  lens.  Second,  since  F,j  is  a  2  x  3  matrix,  there  is  no  unique  inverse  solution  of 
(4)  or  (5)  for  Ax  in  terms  of  AX.  These  difficulties  exist  whether  the  imaging  is  distorted  or  not, 
and  they  are  simply  a  consequence  of  the  effect  of  perspective  (Adrian  1991).  In  deep-field  PTV 
(i.e.  photogrammetry),  the  solution  is  to  use  two  (or  possibly  more)  cameras,  each  giving 
measurements  of  image  locations  in  two  directions  (Nishino  et  al.  1989,  Maas  et  al.  1993,  and 
Guezennec  et  al.  1994).  Then,  the  number  of  measured  image  displacements  exceeds  the  number 
of  unknowns  and  both  x  and  Ax  can  be  found. 

In  light  sheet  PIV,  the  situation  is  a  little  different  because  the  thickness  of  the  light  sheet, 
Azo,  is  normally  very  small  compared  to  the  object  distance,  do,  and  the  width  of  the  field  of  view 
in  the  fluid,  /,  may,  in  certain  cases,  also  be  small  compared  to  dg.  These  conditions  allow  certain 
approximations  to  be  made.  Specifically,  consider  a  PIV  measurement  of  displacement  that  is 
made  in  an  interrogation  spot  centered  at  X/  on  the  image  plane.  It  is  known  (Adrian  1988)  that  a 
PIV  measurement  by  correlation  yields,  on  average,  a  measurement  of  the  velocity  field  of  the 
fluid  averaged  over  a  volume  consisting  of  the  intersection  of  the  interrogation  spot,  projected 
back  into  the  fluid  with  the  thickness  of  the  light  sheet,  Azo  (c.f.  Fig.  1).  This  volume  is  called 
the  interrogation  cell.  The  center  of  the  interrogation  spot  X/  maps  back  into  a  point  that 
intersects  the  mid-plane  of  the  light  sheet  at 


Tt,  ={x,,y,,0),  (8) 

where 


X,  =F(x,). 


(9) 
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This  is  the  center  of  the  interrogation  cell,  and  since  X/  is  two-dimensional  {zi  =  0)  it  is  possible  to 
associate  each  X/  with  a  unique  X/,  thereby  resolving  the  first  difficulty  cited  above. 

The  first-order  relationship  between  an  image  plane  displacement,  AX(X/),  and  object 
plane  displacement,  Ax(x/),  is  derived  by  Taylor  series  expansion  of  Eq.  (4)  about  the  point  X/, 
followed  by  volume  averaging  the  equation  over  the  interrogation  cell,  with  suitable 
transformations  between  x  and  X.  Denoting  averages  over  the  interrogation  cells  in  x  and  X  by 
overbars,  and  neglecting  any  variation  of  the  magnification  matrix  over  the  dimensions  of  the 
interrogation  cell  (justified  when  the  dimensions  of  the  interrogation  cell  are  much  smaller  than 
the  object  distance),  one  finds  the  following  set  of  equations: 

AXi  =  Fij{x,)Axj,  /  =  l,2and;  =  1,2,3,  (10a) 


or 


AX2 

Ax, 

'ax.' 

^2.2 

(10b) 


This  is  the  final  operative  system  of  equations  for  interpreting  PIV  measurements  in  the  presence 
of  distortion.  If  the  distortion  is  severely  nonlinear,  or  if  very  high  accuracy  is  needed  in  the 

relationship  between  AX  and  Ax,  and  justified  by  the  accuracy  of  the  measurements  of  AX ,  this 
relationship  can  be  improved  by  retaining  higher  order  terms  in  the  Taylor  series  expansion.  The 
work  presented  here  assumes  that  a  linear  expansion  is  adequate,  as  is  standard  practice  in  all 
PIV  work  done  to  date. 

For  the  case  of  a  single  camera,  the  system  of  equations  in  (10)  cannot  be  solved  for 
Ax(x/)  as  it  stands  because  there  are  still  only  two  equations  for  three  unknowns.  Common 
practice  in  PIV  with  one  camera  is  to  ignore  the  effect  of  perspective,  as  embodied  in  the  terms 
Fi,3  and  F2,3,  thereby  reducing  Eq.  (10b)  to 


(axJ 


F,  2  Y  Ax, 

Ax, 


2.2 


(11) 


This  is  justified  if  Fi,3  and  F2,3  are  small.  Some  guidance  concerning  the  magnitude  of  the 
perspective  effect  can  be  gained  by  considering  the  expression  for  the  magnification  matrix  at  x/ 
for  perfectly  linear  imaging: 


^,y(X/) 


-d: 


1  0 
0  1 


(12) 


It  is  clear  that  the  perspective  effect  is  small  if  x/  and  yi  are  both  small  compared  to  do  for  all 
values  of  x/  and  y/,  i.e.  over  the  full  field-of-view.  This  condition  is  referred  to  as  paraxial 
imaging  (Adrian  1991).  It  is  expected  that  paraxial  imaging  will  lead  to  small  magnitude  of  the 
perspective  effect  in  the  case  of  distorted  imaging,  as  well,  but  this  must  be  verified  in  individual 
cases. 
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Unfortunately,  it  is  quite  common  for  the  paraxial  condition  to  be  only  weakly  satisfied. 
For  example,  consider  a  typical  laboratory  configuration  in  which  a  100  mm  wide  field-of-view 
is  imaged  at  an  object  distance  of  500  mm.  The  position  x/  varies  from  -50  mm  to  +50  mm, 
implying  that  F\j  varies  between  -0.1  and  +0.1.  Thus,  without  correcting  for  perspective, 
measurements  of  the  in-plane  velocity  components  are  contaminated  by  as  much  as  10%  of  the 
out-of-plane  velocity.  If  the  flow  is  two-dimensional  in  the  Jiy-plane,  i.e.  w  =  0,  there  is  no  error, 
but  if  the  out-of-plane  component  is  comparable  to  the  in-plane  component,  errors  of  order  10% 
are  incurred.  Stereoscopic  imaging  provides  the  only  unconditional  solution  to  the  problem  of 
compensating  for  perspective,  and  it  will  be  considered  in  the  next  section. 


2.2  Stereoscopic  Camera 

For  a  stereoscopic  system,  two  cameras  are  involved  in  deriving  the  displacement 
equations,  which  requires  a  slight  change  in  notation.  The  displacement  of  a  particle  image  on 
camera  c,  where  c  =  1,  2,  is 

AX‘‘>=F‘^»(x  +  Ax)-F‘">(x).  (13) 


Performing  a  Taylor  series  expansion  of  Eq.  (13)  and  volume  averaging  over  the  interrogation 
cell,  similarly  to  what  was  done  in  the  previous  section,  one  finds  the  first-order  relationship 
between  an  image  plane  displacement  on  camera  c,  AX^'^^X/),  and  an  object  plane  displacement, 
Ax(x/),  to  be 

=  F.[f{x,)A^j,  1  =  1, 2  and  7  =1,2, 3.  (14) 


Eq.  (14)  provides  a  system  of  equations  for  each  camera.  Writing  out  both  sets  of  equations  in 
full  £ind  augmenting  them  yields 


or 


AX,^'^  ' 

f 

^i.i 

17(1) 

^1.2 

17(0  \ 

AX<‘> 

r(l) 

^2,1 

pm 

^2,2 

rO) 

^  2,3 

AX<^^ 

17(2) 

^1.1 

77(2) 

^1,2 

AX‘2> 

r(2) 

pi2) 

^2,2 

P2^] 

Ax, 

Axj 

Ax, 


(15a) 


AXa  =  (VF)^Ax, 


(15b) 


where  the  subscript  A  denotes  the  appropriate  augmented  entity.  This  is  the  final  operative 
system  of  equations  for  a  stereoscopic  camera  for  interpreting  PIV  measurements  in  the  presence 
of  distortion. 

Although  it  appears  that  the  system  in  (15)  is  overdetermined,  this  is  not  necessarily  the 
case.  If  no  errors  are  incurred  when  measuring  the  image  displacements,  two  of  the  equations 
will  be  dependent  simply  due  to  the  nature  of  the  optics.  For  example,  if  a  purely  X3- 
displacement  is  observed  by  a  stereocamera,  the  perceived  X2-displacement  on  both  cameras  will 
be  identical,  but  the  perceived  xi -displacement  will  be  in  opposite  directions  on  each  camera. 
Thus,  the  second  and  fourth  equations  in  (15)  will  be  linearly  dependent.  Several  approaches 
could  then  be  used  to  solve  the  system.  For  example,  the  two  dependent  equations  could  be 
added  together  yielding  a  solvable  system.  However,  random  errors  are  incurred  when 


8 


measuring  the  image  displacements,  so  it  may  not  be  the  case  that  two  of  the  equations  will  be 
exactly  dependent.  For  the  general  case  where  the  redundancy  is  not  obvious,  a  least-squares 
solution  should  be  used  to  minimize  the  squared  Euclidean  norm  of  the  residual. 


£  = 


AX/i  -(VF)^Ax 


(16) 


3.  Procedure  for  Single  Camera 

The  purpose  of  the  calibration  is  to  determine  how  positions  in  the  object  plane  map  to 
positions  in  the  image  plane  taking  into  account  all  aberrations.  In  order  to  do  this  it  is  necessary 
to  acquire  images  of  a  calibration  target  whose  location  in  the  object  plane  is  known.  Design  of 
the  calibration  target  is  the  most  important  part  of  performing  a  successful  calibration.  In 
general,  different  experiments  may  require  different  calibration  targets. 

The  target  must  have  some  sort  of  characteristic  marks  on  its  face.  These  marks  and  their 
corresponding  images  will  define  how  positions  in  the  object  plane  map  to  positions  in  the  image 
plane.  The  marks  are  placed  on  a  precise,  evenly  spaced  Cartesian  grid.  The  size  and  shape  of 
the  target  is  constrained  by  the  measurement  container,  but  it  must  be  large  enough  so  that  the 
characteristic  marks  occupy  the  entire  field  of  view  of  the  camera.  In  order  to  maximize  the 
accuracy  of  the  mapping  function  the  field  of  marks  should  extend  to  the  edges  of  the  viewable 
field. 

In  general,  the  grid  on  the  calibration  target  has  two  parameters:  (1)  the  total  number  of 
points,  and  (2)  the  spacing  between  points  (the  size  of  a  grid  cell).  These  parameters  should  be 
chosen  according  to  the  size  of  the  image  acquisition  area  and  the  amount  of  distortion  in  the 
experiment.  Thus,  the  number  of  grid  points  is  proportional  to  the  size  of  the  image  acquisition 
area,  while  the  size  of  a  grid  cell  is  inversely  proportional  to  the  amount  of  distortion  in  the 
experiment. 

A  sometimes  attractive  alternative  to  a  calibration  grid  is  a  single  target  point  that  can  be 
moved  precisely  to  many  different  locations.  For  example,  Kent  et  al.  (1993)  describe  a  target  on 
the  end  of  an  arm  under  control  of  a  stepper  motor.  This  approach  is  very  useful  for  measuring 
the  mapping  function  over  a  deep  volume,  as  in  photogrammetric  PTV. 

Ideally,  one  would  like  to  find  a  mapping  function  F  such  that  given  any  point  x  in  the 
object  plane,  its  corresponding  location  X  in  the  image  plane  could  be  determined  with  complete 
accuracy.  In  reality,  one  cannot  experimentally  determine  an  F  that  will  yield  zero  error  for  an 

arbitrarily  complex  geometry.  Thus,  an  estimate  for  the  mapping  function,  F,  that  has  an 
analytic  representation  will  be  sought  instead.  Many  possibilities  arise  at  this  point,  but  because 
of  its  simplicity,  a  multi-dimensional  polynomial  is  an  obvious  choice.  There  are  several  ways  to 
determine  the  polynomial.  For  example,  one  could  use  the  method  of  least-squares  to  determine 
an  F  that  approximates  the  data  while  minimizing  the  mean  square  error. 


£  = 


2F(x,)-F(xJ 

(xj 


2 


(17) 


where  the  set  {Xc}  represents  the  collection  of  points  on  the  calibration  target,  and  Nc  is  the  total 
number  of  points  in  the  set.  Alternatively,  an  interpolating  spline  could  be  used  to  determine 
local  polynomials. 
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In  all  respects,  the  least-squares  polynomial  is  the  easiest  form  to  implement  and  yields 
very  good  results  unless  the  distortion  is  extremely  severe.  In  these  instances,  a  spline  would  be 
a  better  choice.  A  least-squares  polynomial  with  cubic  dependence  in  x\  and  xz,  but  quadratic 
dependence  in  X3,  will  be  adopted  as  the  mapping  function  estimate  for  the  remainder  of  this 
paper.  That  is, 

F(x)=ao  +a,JC]  +^2X2  -1-33X3  +2iixf  +  a^x^X2  +  Si^xl  +87X3X3  + 

83X2X3 +39X3 +a,oX,^ +3, jXj^Xj +a,2J:,JC2 +^13^2  >  (^8) 

8,4X1  X3  +  815X1X2X3  +  335X2X3  +  837X3X3  +  833X2X3 , 

where  the  a,  are  vector-valued  coefficients  to  be  determined.  It  is  noted,  however,  that  the 
generalization  to  other  methods  of  interpolation  or  approximation  is  trivial  and  should  be  done  if 
it  is  warranted  by  the  severity  of  the  distortion. 

Choosing  F  to  be  a  cubic  function  of  xi  and  X2  was  a  decision  based  purely  on  qualitative 
information.  By  viewing  a  Cartesian  grid  through  various  types  of  distorting  media  that  might  be 
found  in  a  PIV  or  PTV  experiment,  it  was  observed  that  the  curvature  of  the  grid  lines  did  not 
vary  much.  Thus,  it  was  judged  that  a  cubic  polynomial  could  track  the  distortion  produced  in 
the  most  conunon  cases  reasonably  well.  The  justification  for  selecting  quadratic  dependence  in 
X3  is  due  solely  to  the  mathematical  limitation  that  results  from  performing  the  calibration  on 
only  three  X3-planes.  Generalizing  the  procedure  to  four  or  more  planes  would  clearly  allow 
higher-order  terms  in  X3,  but  for  present  purposes,  three  planes  are  adequate. 

Before  beginning  the  calibration  procedure,  one  must  define  a  coordinate  system  in  the 
object  plane  and  the  image  plane.  The  most  common  way  to  define  the  object-plane  coordinate 
system  is  by  specifying  one  of  the  characteristic  marks  on  the  face  of  the  calibration  target  as  the 
origin.  The  locations  of  all  other  characteristic  marks  are  measured  relative  to  this  point. 
Normally,  the  origin  of  the  image  plane  is  defined  by  pixel  (0,0)  on  a  CCD  array  or  the  lower-left 
comer  of  a  photograph.  One  must  also  make  sure  the  experimental  setup  is  in  the  exact 
configuration  it  will  be  in  when  measurements  are  taken. 

The  first  step  in  the  calibration  procedure  is  to  acquire  a  calibration  image.  This  is  simply 
an  image  of  the  calibration  target  as  viewed  by  the  camera  in  the  current  experimental  setup. 
Before  acquiring  the  image,  the  face  of  the  calibration  target  should  be  aligned  as  precisely  as 
possible  with  the  center  of  the  light  sheet.  Once  aligned,  the  calibration  image  is  acquired  and 
stored. 

The  next  step  is  to  determine  the  locations  of  the  characteristic  marks  on  the  calibration 
image  with  respect  to  the  origin  of  the  image  plane.  An  accurate  way  to  locate  the  marks  is 
through  a  template-matching  technique.  This  process  is  performed  in  a  manner  similar  to  the 
correlation  procedure  used  in  PIV.  A  template  image  is  created  that  is  similar  to  the  shape  of  the 
characteristic  marks  on  the  calibration  target.  The  template  and  calibration  image  are  then 
correlated  using  two-frame  cross-correlation,  which  is  followed  by  a  peak-searching  routine  to 
determine  the  pixel  locations  of  the  marks  in  the  correlation  plane.  Sub-pixel  accuracy  can  be 
attained  by  fitting  a  curve,  such  as  a  Gaussian  or  a  quadratic  polynomial,  to  the  peak. 
Alternatively,  for  marks  of  an  appropriate  shape  (e.g.  circular),  one  can  use  a  local  centroid 
calculation,  rather  than  template  matching,  to  determine  the  locations  of  the  marks.  Once  the 
locations  of  the  target  points  in  the  image  and  object  planes  have  been  obtained,  the  coefficients 
of  the  mapping  function  in  Eq.  (18)  can  be  determined  using  a  least-squares  procedure.  Note  that 
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for  a  single  camera  the  calibration  image  is  acquired  at  only  a-single  JC3  location  (presumably  at  xj 
=  0).  In  this  case,  all  terms  in  Eq.  (18)  involving  ;C3  are  zero,  and  the  corresponding  coefficients 
are  not  determined  by  the  calibration.  This  is  acceptable  because,  as  noted  above,  the  single 
camera  equations  must  ignore  out-of-plane  motion  anyway. 

Now  that  the  system  has  been  calibrated,  one  may  proceed  normally  and  acquire  images 
using  the  conventional  PIV  technique.  Once  recorded,  each  image  is  then  interrogated  on  a  set  of 
points  {X/}  in  the  image  plane  that  is  determined  by  the  interrogation  software.  The  image  plane 
displacements  are  then  validated  and  filtered  using  appropriate  software  tools.  Note  that  many 
interrogation  software  systems  provide  the  user  with  very  limited  control  over  the  grid. 
Commonly,  one  can  manipulate  only  a  few  parameters,  such  as  the  location  of  the  grid  in  the 
image  plane  and  the  number  of  interrogation  points  to  use,  and  specifying  an  arbitrary  grid  for 
interrogation  is  not  usually  possible.  Even  if  it  were,  locations  in  the  image  plane  are  not 
necessarily  meaningful.  One  would  prefer  to  be  able  to  request  particle  displacements  at 
locations  in  the  object  plane  for  which  the  coordinates  of  physically  meaningful  locations  are 
known. 

Because  it  is  known  how  positions  in  the  object  plane  map  to  positions  in  the  image  plane 
from  the  mapping  function  generated  during  the  calibration  procedure,  it  is  possible  to  calculate 
particle  displacements  at  any  location  in  the  fluid  within  the  domain  of  the  mapping  function. 
The  collection  of  points  in  the  object  plane  where  displacements  are  desired  is  denoted  by  {x^} 
and  their  corresponding  images  by  {X^}.  The  image  displacements  previously  determined  on 
{X/}  during  the  interrogation  stage  are  interpolated  onto  {X^}  using  an  appropriate  method,  such 
as  bilinear  interpolation.  Note  that  two  potential  problems  can  arise  during  this  procedure.  First, 
the  image  of  a  particular  may  lie  outside  the  image  recording  domain  of  the  camera,  i.e.  the 
CCD  chip  or  the  film  frame.  Second,  when  interpolating  image  displacements  onto  a  particular 
Xg,  one  or  more  of  its  neighboring  grid  points  may  not  contain  data.  In  both  cases,  the  x^  in 
question  is  discarded.  At  this  point  Eq.  (1 1)  is  solved  for  the  particle  displacements  in  the  fluid. 

Fig.  3  is  a  flow  chart  that  illustrates  the  processes  described  above  to  compute  the  particle 
displacements  in  the  fluid.  Process  A  represents  the  steps  taken  during  the  PIV  experiment  and 
ensuing  interrogation  process.  Each  step  in  Process  A  can  be  made  within  any  of  several 
different  standard  PIV  processing  software  packages,  including  commercial  software.  Process  B 
represents  the  procedure  by  which  image  displacements  computed  during  process  A  are 
interpolated  onto  the  mapping  of  the  user-defined  grid.  Process  C  represents  the  calibration 
procedure.  The  results  of  processes  B  and  C  serve  as  input  to  Eq.  (11),  from  which  the  particle 
displacements  on  the  user-defined  grid  are  calculated.  Note  that  defining  the  grid  {x^}  in  the 
object  plane  and  mapping  it  onto  the  image  plane  circumvents  the  nonlinear  inverse  problem  of 
solving  for  object  points  in  terms  of  image  location. 

4.  Procedure  for  Stereoscopic  Camera 

The  calibration  procedure  for  a  stereoscopic  camera  is  very  similar  to  that  for  a  single 
camera  (c.f.  Fig.  4).  However,  because  it  is  necessary  to  know  how  the  magnification  in  the  x^- 
direction  varies,  one  must  acquire  multiple  calibration  images  at  planes  normal  to  the  light  sheet. 
Since  a  polynomial  is  being  used  to  estimate  the  mapping  function,  the  number  of  images 
acquired  at  different  xj  locations  places  an  upper  limit  on  the  order  of  the  x^  terms.  Keep  in  mind 
that  a  least-squares  method  is  being  used  so  that  one  could  use  more  calibration  images  than  the 
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Figure  3.  Flow  chart  for  a  single  camera.  Process  A  represents  the  standard  PIV  acquisition  and 
interrogation  procedure.  Process  B  represents  the  procedure  for  computing  displacements  in  the  image 
plane  corresponding  to  user-defined  locations  in  the  fluid.  Process  C  represents  the  calibration  procedure. 
The  results  of  processes  B  and  C  are  used  to  compute  the  particle  displacements  at  user-defined  locations 
in  the  fluid. 
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CALIBRATION  TARGET 


Figure  4.  Schematic  of  the  calibration  procedure  for  a  stereoscopic  camera.  The  mapping  function  for 
each  camera  is  valid  over  the  volume  of  size  dx  x  dy  x  dz  traced  out  by  the  calibration  target  during  the 
calibration  procedure. 
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desired  order  of  the  xa-terms  in  order  to  reduce  the  mean  square  error  of  the  estimate.  At  a 
minimum,  the  second  order  polynomial  in  xt,  requires  calibration  data  at  three  X3-locations 
separated  by  a  distance  comparable  to  the  thickness  of  the  light  sheet. 

Once  the  calibration  images  for  all  planes  have  been  acquired  for  each  camera,  the 
locations  of  the  characteristic  marks  in  each  image  are  determined  with  respect  to  the  origin  of 
the  corresponding  image  plane.  The  same  procedure  described  in  the  previous  section  for 
locating  the  images  of  the  calibration  marks  is  applied  to  each  image.  Enough  information  is 
now  available  to  determine  the  coefficients  in  Eq.  (18)  twice  —  once  for  each  camera.  As  shown 
in  §2.2,  the  resulting  two  mapping  functions  are  all  that  is  necessary  to  recombine  a  pair  of  2-D 
velocity  fields  into  a  single  planar  3-D  velocity  field. 

Once  the  system  has  been  calibrated,  one  can  proceed  normally  and  acquire  images  using 
the  conventional  PIV  technique.  The  images  for  each  camera  are  interrogated  on  a  set  of  points 

in  the  corresponding  image  plane  that  is  determined  by  the  interrogation  software.  The 
limitations  mentioned  in  the  previous  section  concerning  user-control  of  the  grid  apply  here,  as 
well.  The  ability  to  specify  a  user-defined  grid  in  the  fluid,  as  in  the  single  camera  case,  is 
applicable  to  the  stereoscopic  camera,  also.  The  two  problems  that  may  arise  during  the 
interpolation  of  image  displacements  onto  can  now  occur  for  each  camera  individually.  If 

a  particular  \g  has  been  deleted  because  of  a  problem  that  occurred  during  interpolation  on  either 
camera,  stereo  recombination  cannot  be  performed  at  this  location  and  a  hole  will  appear  in  the 
final  3-D  velocity  field. 

Fig.  5  is  a  flow  chart  that  illustrates  the  processes  described  above  to  compute  3-D 
particle  displacements  in  the  fluid.  The  processes  outlined  in  Fig.  3  are  performed  for  each 
camera.  The  results  of  processes  Bl,  B2,  Cl,  and  C2  are  used  to  form  the  augmented  system  in 
Eq.  (15),  and  the  3-D  particle  displacements  on  the  user-defined  grid  are  calculated.  It  should  be 
noted  that  unlike  the  stereoscopic  PIV  recombination  algorithm  given  in  Prasad  and  Adrian 
(1993)  in  which  the  image  plane  interrogation  grid  and  corresponding  image  displacements  were 
inversely  mapped  into  the  fluid  followed  by  the  determination  of  a  common  Cartesian  grid  to 
interpolate  displacements  onto,  the  algorithm  presented  here  requires  no  such  procedure.  The 
entire  process  is  simplified  by  mapping  the  user-defined  grid  {x^}  onto  each  image  plane  and 
performing  the  interpolation  of  displacements  therein. 

5.  Application  and  Evaluation  of  Calibration  Procedure 

In  this  section,  we  present  the  results  of  applying  the  calibration  procedure  to  both  a 
single  camera  and  a  stereoscopic  camera.  In  each  case,  a  computer-generated  image  was 
displaced  by  a  known  amount  in  the  object  plane  and  the  corresponding  image  displacements 
were  computed  using  TSI  INSIGHT  PIV  software.  Mapping  functions  were  used  to  compute  the 
fluid  displacement  from  this  information  and  errors  are  plotted.  Unless  otherwise  mentioned,  the 
uncertainty  in  any  experimental  displacements  of  the  image  is  ±  5  pm. 

5.1  Single  Camera 

The  setup  for  the  single  camera  experiment  is  shown  schematically  in  Fig.  6.  In  order  to 
test  the  robustness  of  the  calibration  procedure,  a  distorting  medium  was  purposely  placed 
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Figure  5.  Flow  chart  for  stereoscopic  configuration.  Processes  A,  B,  and  C  from  Fig.  3  are  performed 
for  both  cameras.  An  augmented  system  is  formed  and  solved  for  the  3D  particle  displacements  at  user- 
defined  locations  in  the  fluid.  Note  that  (VF)a*  is  the  generalized  inverse  of  (VF).^. 


Distorting 

Medium 


Figure  6. 


between  the  lens  and  the  object  plane.  Two  lenses  of  focal  length  450  mm  placed  in  series  were 
used  to  provide  the  distortion.  Fig.  7(a)  illustrates  the  effect  of  the  distortion  on  a  uniform 
Cartesian  grid  placed  in  the  object  plane. 

The  calibration  target  used  in  this  experiment  was  a  uniform  12  x  12  Cartesian  grid  of 
white  circles  on  a  black  background.  The  diameter  of  each  circle  was  1  mm,  while  the  spacing  of 
the  circles  was  5  mm  in  x  and  y.  The  calibration  procedure  was  performed  and  the  mapping 
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function  was  generated.  The  RMS  error  in  the  mapping  function  was  calculated  to  be  0.53  pixels 
(corresponding  to  an  average  RMS  error  of  6.9  x  10”"  mm  in  the  object  plane).  It  should  be 
noted  that  the  RMS  error  is  influenced  by  not  only  the  ability  of  the  least-squares  polynomial  to 
fit  the  data,  but  also  by  the  accuracy  with  which  image  locations  of  the  characteristic  marks  are 
measured  during  calibration.  Prasad  et  al.  (1992)  have  shown  that  the  random  error  associated 
with  locating  particle  image  centroids  in  PIV  is  approximately  7%  of  the  particle  image  diameter. 
The  random  error  associated  with  the  least-squares  fit  is  approximately  6%  of  the  particle  image 
diameter.  Thus,  the  accuracy  of  the  curve  fit  is  comparable  to  the  accuracy  for  locating  particle 
image  centroids,  from  which  it  is  inferred  that  the  third-order  polynomial  was  adequate  to 
represent  the  distortion  in  this  experiment. 

To  test  the  algorithm  for  computing  accurate  fluid  displacements,  a  computer-generated 
pattern  of  randomly  located  dots  was  placed  on  a  translation  stage,  and  an  image  was  acquired  at 
both  a  reference  position  and  after  a  displacement  of  Ax  =  (0.500,  1.000,  0.000)  mm.  The  two 
images  were  then  interrogated  using  two-frame  cross-correlation  to  determine  the  displacements 
of  the  particle  images,  and  the  resulting  vector  field  was  validated  and  filtered. 


Figure  7.  (a)  Distorted  image  of  a  uniform  Cartesian  grid,  and  (b)  undistorted  image  obtained  using  the 
mapping  function  from  the  calibration.  The  region  of  (a)  that  has  been  undistorted  is  indicated  by  the 
dashed  rectangle. 
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Figs.  8  and  9  are  contour  plots  of  the  error  in  the  computed  fluid  displacement  obtained 
using  the  displacement  equations  for  a  single  camera.  The  error  has  been  nondimensionalized  by 
the  size  of  the  mapping  function  domain  in  the  appropriate  direction.  In  this  experiment,  the  size 
of  the  domain  was  55.8  mm  in  both  x  and  y.  This  nondimensionalization  was  chosen  in  order  to 
show  how  much  the  calculated  fluid  displacements  are  in  error  relative  to  the  size  of  the  region 
being  mapped.  The  largest  error  in  both  figures  is  approximately  1/1000  the  size  of  the  mapping 
function  domain.  It  does  not  appear  that  any  trends  are  present  in  the  error,  implying  that  random 
error  is  the  dominating  source.  However,  there  are  several  areas  in  both  figures  where  the  error 
appears  to  be  correlated. 

As  noted  in  §2,  it  is  possible  to  compute  an  inverse  mapping  function  for  the  single 
camera  since  each  X/  is  associated  with  a  unique  X/.  Using  the  same  calibration  image  previously 
acquired,  the  inverse  mapping  function  for  this  setup  was  generated,  and  its  RMS  error  was 
calculated  to  be  1.5  x  10"^  mm.  The  inverse  mapping  function  was  used  to  “dewarp”  the  grid 
shown  in  Fig.  7(a),  and  the  dewarped  grid  appears  in  Fig.  7(b). 


2.5x1 0'^ 
1.5x10-^ 
5.0x1  O'* 
-5.0x10' 
-1.5x10'' 


Figure  8.  Error  between  the  true  jr-displacement  in  the  object  plane  and  the  computed  value.  The  error 
has  been  nondimensionalized  by  55.8  mm,  which  is  the  size  of  the  mapping  function  domain  in  the  x- 
direction. 


17 


-20  -10  0  10  20 

x(mm) 


3.5x10'®^ 
2.5x10'°^ 
1.5x10"^ 
5.0x1 0"“ 
-5.0x10-“ 
-1.5x10-“ 


Figure  9.  Error  between  the  true  ^-displacement  in  the  object  plane  and  the  computed  value.  The  error 
has  been  nondimensionalized  by  55.8  mm,  which  is  the  size  of  the  mapping  function  domain  in  the  y- 
direction. 


5.2  Stereoscopic  Camera 

The  stereoscopic  camera  is  shown  schematically  in  Fig.  10.  For  a  detailed  description  of 
the  camera  configuration,  see  Liu  et  al.  (1997).  In  order  to  test  the  robustness  of  the  mapping 
function  approach  for  stereo  recombination,  a  distorting  medium  was  purposely  placed  between 
the  lenses  and  the  object  plane.  A  water-filled  fish  tank,  approximately  305  mm  deep,  was  used 
to  provide  the  distortion. 

The  calibration  target  used  in  this  experiment  was  a  rectangular  block  of  aluminum  with 
81  circular  holes  drilled  into  its  face.  These  holes  were  arranged  on  a  uniform  9x9  Cartesian 
grid.  The  diameter  of  each  hole  was  0.5  mm,  while  the  spacing  of  the  holes  was  9  mm  in  x  and  y. 
A  light  source  placed  behind  the  target  was  used  to  illuminate  the  holes.  The  calibration 
procedure  was  performed  and  the  mapping  function  for  each  camera  was  generated.  During  the 
calibration,  the  target  was  placed  at  xj  locations  of  0.000,  0.500,  and  -0.500  mm.  The  RMS 
errors  in  the  mapping  functions  for  cameras  1  and  2  were  1.1  pixels  and  1.2  pixels,  respectively. 
These  correspond  to  an  average  RMS  error  of  4.5  x  10“^  mm  and  5.1  x  10“^  nun  in  the  object 
plane,  respectively. 
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Figure  10.  Schematic  of  the  stereoscopic  camera. 

To  test  the  recombination  algorithm,  a  computer-generated  pattern  of  randomly  located 
dots  was  placed  on  a  translation  stage,  and  an  image  was  acquired  at  a  reference  position  and 
after  a  displacement  of  Ax  =  (0.000,  1.000,  -0.400)  mm.  The  two  images  for  each  camera  were 
then  interrogated  using  two-frame  cross-correlation  to  determine  the  displacements  of  the  particle 
Tmages.  The  resulting  2-D  vector  fields  were  then  subjected  to  vector  validation  and  filtering. 
Finally,  the  stereoscopic  recombination  algorithm  was  applied  and  the  planar  3-D  velocity  field 
was  generated. 

Figs.  11-13  are  contour  plots  of  the  error  in  the  computed  fluid  displacement  obtained 
using  the  displacement  equations  for  a  stereoscopic  camera.  The  error  has  been 
nondimensionalized  by  the  size  of  the  mapping  function  domain  in  the  appropriate  direction.  In 
this  experiment,  the  size  of  the  domain  was  72  mm  in  both  x  and  y  and  1  mm  in  z.  As  for  the 
single  camera,  it  appears  that  no  trends  are  present  in  the  error,  but  there  are  still  correlated 
regions. 

6.  Summary 

A  procedure  for  calibrating  both  conventional  and  stereoscopic  PIV  systems  has  been 
presented.  The  results  of  the  calibration  yield  one  or  more  mapping  functions  that  provide  a 
convenient  and  accurate  way  to  correct  for  distortion  in  the  experiment  when  determining 
particle  displacements.  In  particular,  it  is  possible  to  accurately  determine  the  magnification 
matrix  that  relates  image  plane  displacements  to  object  space  displacements.  In  addition,  a 
stereoscopic  camera  is  automatically  registered  during  the  calibration  procedure,  which  frees  the 
user  from  the  task  of  mechanically  registering  the  camera.  This  permits  rapid  and  easy  setup  of 
angular-offset  stereoscopic  cameras.  The  analysis  provides  accurate  equations  by  which  one 
solves  for  the  2-D  or  3-D  displacement  of  the  fluid  and  it  permits  the  use  of  arbitrary  user- 
defined  grids  in  the  fluid. 
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Figure  11.  Error  between  the  true  x-displacement  in  the  object  plane  and  the  computed  value.  The  error 
has  been  nondimensionalized  by  72  mm,  which  is  the  size  of  the  mapping  function  domain  in  the  x- 
direction. 
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Figure  12.  Error  between  the  true  y-displacement  in  the  object  plane  and  the  computed  value.  The  error 
has  been  nondimensionalized  by  72  nun,  which  is  the  size  of  the  mapping  function  domain  in  the  y- 
direction. 
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Figure  13.  Error  between  the  true  z-displacement  in  the  object  plane  and  the  computed  value.  The  error 
has  been  nondimensionalized  by  1  mm,  which  is  the  size  of  the  mapping  function  domain  in  the  z- 
direction. 
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(217)  244-5707  -  Pax 


Quote  No.  R2-189-2  (Page  1  of  2) 

Quote  Date  January  31,  1995 
Customer  Ref. 

Due  Date 

Terms:  Net  30  Days 

Firm  Quote  for _ 90 _ Days 

Prices  are  F.Q.B.  St.  Paul,  MN  U.S.A. 


Ul"™  Quantity 


TSI 

Model  No. 
600031 


600041 


Description 


Insight  PC  Window  based  PIV  Analysis  Software 
Acquisition  and  Processing  Module  with 

*  Support  for  both  film  based  and  video 
based  PIV  data 

*  Selectable  pixel  resolution  and  user 
configurable  grid  size 

*  2-D  auto  and  cross  correlation 

*  Real-time  control  and  manipulation  on  image 
file 

*  Open  Image  format 

*  1-5  vectors /sec  at  486  based  computer 

*  10-40  vectors/sec  with  array  processor 
Post  Processing  Module  with 

*  Automated  Data  Validation  and  Recovery 

*  Multiple  Window  Display 

*  Velocity  Averages  and  Derivatives 

*  Printer  Support 

Insight  Unix  Workstation  based  PIV  Analysis 
Software  Acquisition  and  Processing  Module  wit 

*  Support  for  both  film  based  and  video  based 
PIV  data 

*  Selectable  pixel  resolution  and  user 
configurable  grid  size 

*  2-D  auto  and  cross  correlation 

*  Real-time  control  and  manipulation  on  image 
file 

*  Open  Image  format 

*  10-100  vectors/sec  output 
Post  Processing  Module  with 

*  Automated  Data  Validation  and  Recovery 

*  Multiple  Window  Display 

*  Velocity  Averages  and  Derivatives 

*  Printer  Support 


Unit  Extended 

Price  Total 


$  8,000  $  8,000.00 


10,000 


10,000.00 


Delivery  Terms 


6  weeks  ARO 


This  Quotation  is  subtect  to  the  warranties,  disdaimers  and  all  other  terms  and  cortditions  sat  forth  on  the 
reverse  side  hereof  and  to  no  others.  Seller  reserves  the  right  to  chartge  prices  effective  on  any  new  orders, 
provided  Seller  notifies  in  writing  those  with  cufrendy  valid  Quoutions  prior  to  any  order  being  placed.  This 

Accepted  By  Buyer 

By: - - - 


quotation  shall  become  an  agreement  binding  upon  the  Buyer  and  the  Seller  when  accepted  by  the  Buy 
and  subsequently  accepted  by  an  authonzed  representative  of  the  Seller  at  the  Seller’s  home  office  ar 
thereupon  shaU  constitute  the  entire  agreement  between  the  parties. 


Martin  Abbott  ^  V  | 

Sales  Manager  -  Midwest  Region 


Form  No.  Q100A 


Quotation 


TSf.l 


rSI  incorponted 

Fluid  Mechanics  Instalment  Division 
500  Cardigan  Road  612  490  281 1 

P.  O.  Box  64204  Telex  6879024 

St.  Paul,  MN  55164  USA  Fax  612  490  3824 


